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Abstract 

In the present work, a new time-dependent exchange theory is presented wherein the sym- 
metry constraints, on a multi-electron wavefunction, are properly accounted for. In so 
doing, the equations of motion, incorporating the required symmetry, are derived and a 
solution algorithm employing an implicit split-operator procedure is described. A tech- 
nique (using an orthonormalization transformation and a unitary rotation), for explicitly 
enforcing the required constraints, which render the computations tractible and provide 
for non-redundant time evolution, is also presented. This amounts to the calculation of the 
appropriate numerically determined guage. The invariance of the derived orbital equations 
of motion with respect to the transformations is explicitly demonstrated. 
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I. INTRODUCTION 



The time-dependent Schrodinger equation (TDSE) 

i^^^H^ (1) 
dt ^ ' 

describes the dynamics of quantum mechanical systems, and in particular is applicable 
to systems consisting of A^-electrons. The present application is restricted to atomic and 
molecular systems with all nuclei fixed in space (assumed infinitely massive). Note that 
atomic units are used throughout this work.[l] In these units, % = rUe = e = 1, where h is 
Planck's constant divided by 27r, rUe is the electron mass, and e is the electron charge. ^ is 
the system wavefunction and is a function of time (t) and of the coordinates of the particles 
making up the system. H is the Hamiltonian operator and consists of the kinetic energy 
operators of all of the particles in the system, plus the interaction potential between the 
particles as well as any external potential. In general, H might explicitly depend on time 
through the external interaction potential, and while this would present no fundamental 
complication, the present work assumes H does not explicitly depend on time. Also, in 
order to make the essential points of the present work, the number of electrons will be 
restricted to two. This would seem at first to be a drastic reduction in complexity, which 
of course it is, but the essential points can be made most clearly for two electron systems, 
and indeed, two electron systems, such as the hydrogen molecule {H2)., electron-hydrogen 
atom (e + H) scattering, and the helium atom (He), are important systems. 

The TDSE represents an initial value problem such that if the value of the wavefunc- 
tion at t = is specified, and if the TDSE can be accurately solved, then the quantum 
mechanical dynamics will be encoded in the solution (wavefunction) at t = oo or an approx- 
imation thereunto. In section II, the multiconfigurational time-dependent exchange theory 
will be presented in a form which is specialized to two-electron systems. In addition, since 
much of the present development is herein presented for the first time, only one orbital per 
particle will explicitly be considered. It should be noted that the present formulation can 
be viewed as a modification of the multiconfigurational time-dependent Hartree (MCTDH) 
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theory of Manthe, Meyer, and Cederbaum,[2,3] appropriate for fermions, and as such, it 
may be generalized to several fermions and multi-function representations of each one. But 
in addition, the present paper presents an explicit prescription for enforcing the required 
constraints on the time-dependent orbitals such that "non-rcdundant time evolution" is 
guaranteed-this appears to be a new contribution, along with the new explicit treatment 
of the exchange symmetry and the very compact form of the EOM. Section II. A derives 
the equations of motion (EOM) assuming orbital orthonormality and non-redundant time 
evolution; Section II. B presents the derivation of the transformations-(l) an ortho normal- 
ization matrix and (2) a unitary rotation matrix, which is the solution of a first order 
differential equation in time. The application of these two orbital transformations results 
in orbital orthonormality and non-reduntant orbital time evolution; Section II. C demon- 
strates the invariance of the EOM with respect to the two transformations; Section III 
describes a solution algorithm using an implicit split-operator procedure (ISOP)[4,5]. Fi- 
nally, section IV will present the conclusions. 



II. TIME-DEPENDENT ELECTRON EXCHANGE THEORY 

II. A Equations of Motion 

As specified above, the present development is applied to two-electron systems and 
thus the TDSE is represented as 

z^*(l,2) = ^(l,2)*(l,2) (2) 

where the notation for ^ stands for ^^(1,2) = '^{fi,f2',t), and H{1,2) = H{fi,'r2) such 
that the time dependence is assumed for all wavefunctions and orbitals (to be defined 
below) -the Hamiltonian however is assumed to have no explicit time dependence. For a 
two-electron system, the wavefunction factors into a space part times a spin part and the 
total wavefunction must be completely antisymmetric with respect to electron exchange. [6] 
The spin states are either singlets or triplets. The Hamiltonian is spin independent at the 
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level of theory under consideration. More explicitly, the Hamiltonian is given by 



H{1, 2) = ^o(l) + ^o(2) + V{1, 2) = Ho{l, 2) + V{1, 2) (3) 

where Ho{j) = fo{j) + Vo{j) and fo{j) = -^^f.- Note that if V{1, 2) is symmetric, then 
i^(l,2) is also. 

For a two-electron system, the spatial wavefunction, corresponding to a spin singlet 
(^^(+))^ is symmetric with respect to exchange, while the spatial part is antisymmetric 
for a spin triplet (^f^")). The spatial wavefunction then depends on the spatial parts of 
the two electrons and on time. It resides in a six dimensional (6D) coordinate space. If 
it were practical, the full six spatial dimensional wavefunction would be propagated in 
time. However, this is not practical at the present stage of computer technology-it may 
be in the near future however. Certainly, for systems composed of three or more electrons, 
the full dimensional solution is not available and will not be for the forseeable future. 
Thus, a decomposition into a direct product of three dimensional (SD) subspaces (one 
for each electron) is the advisable procedure. Then, the spatial part of the two-electron 
wavefunction may then be expanded, at the minimal expansion length, in the manner of 
[2,3], as 
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^ U (4) 

The A's are purely time-dependent coefficients multiplying each two-electron configuration 
function. The superscript {(f), ±) indicates that the one-electron orbitals, labeled by (j), are 
used, and the permutation symmetry is singlet or triplet (— ). Thus, the singlet spatial 
wavefunction is a superposition of three two-electron configurations ^^'^''^^ and the triplet 
has one two-electron configuration function ^[^'~\ The minimal length of these expansions 
is determined by the required invariance of ^^^^ with respect to arbitrary rotations among 
the one-electron orbitals[2,3] comprising the $5'^''^^'s (see below). ^^■'^'^^ : j = 1,2,3 is 
symmetric with respect to electron exchange and $('^'~) is antisymmetric with respect to 
electron exchange. There are two one-electron orbitals required to describe the two-electron 



configurations: they are labeled : j = 1,2. 

$('^'+)(l,2) (2) 



(1,2)=^ (1)4+) (2) + (2)4+) (1)] 
4^'+) (1,2) =4+) (1)4^) (2) 

*J^'"^(i,2)=^[4"^(i)4"^(2)-4"^(2)4"^(i)] 

Consistent with the decomposition of the 6D space into a direct product of 3D sub- 
spaces, the objective is to derive two coupled time-dependent 3-D equations for the one- 
electron orbitals. In addition, time-dependent equations for the purely time-dependent 
coefficients need to be derived (three for the singlet and one for the triplet-see Eq. (4)). 

Now in general, the 0's are not necessarily computationally orthonormal because of 
numerical inaccuracies. If Eqs. (5) are substituted into Eqs. (4), which are then substituted 
into Eq. (2), aside from nonzero off-diagonal, and non-unit diagonal, overlaps 5"^^^^^'+) = 
< 4+'*|4^'' >) a set of nonzero "derivative overlaps" given by D[p^^ =< (f)\'^^ \4>j^^ >, 
where the 'over-dot' represents a time derivative, will appear in the equations. One of 
the principal distinctions of the MCTDH theory of [2,3] is the use of the purely time- 
dependent coefficients, as indicated in Eqs. (4), as contrasted with, for example, a time- 
dependent Hartree-Fock theory (TDHF) [7], wherein the time-dependence is described 
solely by the one-electron orbitals. As pointed out in [2,3], this additional set of time- 
dependent coeficients produces a redundant description which allows for the incorporation 
of certain constraints. To see this, consider a new set of one-electron orbitals, '4'j^\ related 
to the old set by 

Assuming b is unitary and purely time-dependent, (bb''' = h^h = 1), this can be written 

~_(±) ~J±) 

v;(±) =b(±)0(±) 



where f indicates hermitian conjugate. Assuming Eq. (7) can be inverted (where the tilde 
indicates vector or matrix transpose) as per 

* =^ (8) 

and inserted into Eqs. (5) and then into Eqs. (4), it will be see that Eqs. (4) can be 
written as 



where 
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*W(1,2) = J] •+)$J^'+)(1,2) 

j=i 



$i^'+)(l,2) =4+\l)4+\2) 

^^^'^^ [4-"^ (1)4^^ (2) + 4-"^ (2)4-"^ (1)] 

$^^'+)(l,2) =4+\l)4-'\2) 

$(^'-)(i,2) =^[4-\i)4-\2) - 4-\2)4-\i)] 



(9) 



(10) 



and 



^1 —^1 L'^ii ^^22 '^21 '^12 J 

The significant point is that the four purely time-dependent matrix elements defining b*^^^ 
are completely arbitrary because of the invariance of Eqs. (9) and (10) with respect to the 
transformation described by Eq. (7). 

The utility of the arbitrariness of b*^^) lies in its use to fix four matrix elements 
involving the i/j^s. As a rationale for this, consider that there are four overlaps S^^'"^^ = 
< 4"^^ 14"^^ ^ ^^^^ derivative overlaps D^p"^^ = < 4'^^\4^^ Thus there are eight 
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such matrix elements. However, only four are independent. If the four choices are made 
corresponding to i:>Jf '^^ = 0, i:>^t'^^ = 0, i:>5t'^^ = 0, and S^f'^^ = 0, then it can be easily 
shown that the other four matrix elements are fixed. Thus, I^g'^^ = 0, = d^^j for 

all Actually, the diagonal overlaps are arbitrary constants which may be set to one. 
In section II. B below, a systematic way of implementing these constraints will be given. 
In so doing, it will be seen that the enforcement of orthonormality and non-reduntant 
time-evolution, must be done in two different transformation steps, separated in sequence 
by the time propagation over the time interval At. 

The derivation of the EOM proceeds by a projection method. A critical property in the 
derivation by projection of the EOM for the ^'s is the orthonormality of the two-electron 
configuration functions such that 

«$(V"±)|$(.^'±)»=5.., (12) 

where the double brackets << | >> represents 6D integration over the coordinates of both 
particles, and the zero values of the derivative overlap matrix elements of the two-electron 
configuration functions 

«$f'±)|^(.^'±)»=0. (13) 

These relations follow immediately from the orthonormality of the ip^s and the zero values 
of the one-electron derivative overlap matrix elements. The procedure is thus to substitute 
Eqs. (10) into Eqs. (9), and then to use the result in Eq. (2). Then project from the left 
by $^ and integrate over the coordinates of particles one and two. The result is 

Ar(±) 

= -z ^ « '±) » (14) 

where N^"^^ — 3/1. for the singlet case (-1- sign), j, : 1, 2, 3. For the triplet case, (— sign), 
j,k:l. 

In order to derive the EOM for the time-dependent one-electron orbitals, the singlet 
and triplet cases are considered separately. To derive EOM for the one-electron orbitals, a 

6 



similar projection (as above for the $'s) is empfoyed, except now the one-electron orbitals 
are used. First, project from the left by V'i^'*(l)* on Eq. (2) and integrate over the 
coordinates of particle one. Then use the overlap and derivative overlap constraints (in 
the 'i/;-basis) to obtain 

A^'^W {2)+^A^^^^W (2) + (2) + (2) 

3 (15) 

where < || >i indicates an integral over the coordinates of particle one. 

Secondly, project from the left by ?/;2^''(2)* on Eq. (2) and integrate over the coordi- 
nates of particle two. Then use the overlap and derivative overlap constraints, again, to 
obtain 

-L4^'+)^(+)(i)+4^'+)4+)(i) + i=4^'+Vi^^(i) +i?'^V^^^(i) 

3 (16) 

i=i 

where < 1 1 >2 indicates an integral over the coordinates of particle two. 

Exactly the same two projections arc done for the triplet case as was done for the 
singlet case. The results are (in analogy to Eq. (15)) 



^A^^'-^4-\2) + i=i(^'-)v,(-)(2) = -i< 4-\l)\H\^?'-\l,2) >i A^^'-^ (17) 
and (in analogy to Eq. (16)) 



-^A[^'-m-\l) + i=4^'-Vi-^(l) = -i< 4"^(2)|^|$S^'-)(1, 2) >, Ai'^'-) (18) 

The three EOM are thus given by Eqs. (14), (15), and (16) for the singlet case, and by 
Eqs. (14), (17), and (18) for the triplet case. 



Plugging Eqs. (10) into Eqs. (15) and (16), then plugging Eqs. (10) into Eqs. (17) 
and (18), and then reversing the coordinate labels in Eq. (15) and Eq. (17) (e.g. 1 ^ 2), 
and finally writing in matrix form, results in 



* (1) = -Z h(^'±)(l) *^^'±)(1) (19) 

where the symbols are defined differently for the singlet and triplet cases, and where 



*(^'±)(1) = a('^'±) V^W(l) (20) 
Now, the h and a-matrices are defined as 

h(^'±)(l)= , (21) 



and, for the singlet 



and for the triplet 



75^2 ^3 



Finally, the h^]^-"^^ matrix elements are defined by 

^(2)1^(1, 2)|^f)(2) >2 (24) 

From here on, the actual orbital coordinate labels will be dropped, except in several cir- 
cumstances where it is convenient to exhibit them for clarity. 

II. B Orbital Transformations 

In the derivation of the EOM [Eqs. (14,19)] using the V'-basis, it has been assumed 
that 8^'^'='=) = 1 and D*^'^'='=^ = 0. However, just because these two assumptions have been 
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made in the deriving the EOM, does not automatically result in the enforcement of the 
constraints expressed in the two assumptions. It is well known that if the multielectron 
time-dependent Schrodinger equation is solved exactly as an A^-electron problem, without 
the orbital direct product ansatz, and if the wave function at t = is orthonormal, then, 
in principle, the wave function should remain orthonormal for all time. However, compu- 
tational errors will inevitably accrue and destroy this orthonormality. For an orbital direct 
product decomposition, this orthonormality assumption is in principal, still valid, subject 
to numerical inaccuracies. If this orbital orthonormality is lost because of numerical in- 
accuracies, it is not sufficient to just re-orthonormalize without appropriately modifying 
the other terms which appear in the EOM-in fact, the EOM should be invariant with 
respect to this orthonormalization transformation. It is one of the salient features of the 
present work, that this invariance is explicit for the EOM derived herein. This invariance 
is demonstrated below. Also, a procedure for enforcement is described. 

It can be easily shown that if in some basis {ip for example) S'^'^'^^ = 1, then the 
D-matrix is anti-hermitian (e.g. D'^^'='=^ -|-D(^'='=^^ = 0). On the other hand, if 0*^^'=*=^ = 0, 
then D(^'='=)^ = and S^^^'^^ = C^'^'^) where C^'^'^) is a constant matrix, not necessarily 
1. As indicated above, however, orthonormality and a null D-matrix are consistent with 
each other. 

II.B-1 Orthonormalization Transformation 

This section is concerned with the following overlap matrices: 

gix, J ^< ^ I ^ > 

A symmetric orthonormalization procedure [8] is utilized, at the R-end of the At interval 
(note that the transformation matrix X^^^ is not unitary, however, note that X^^^ = 
X'^='=^^), to go from the x-set to the 0-set: 

^ = E xi^^^S^ (26) 

V 

or in matrix notation 
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^ ^ (27) 



where 



X(±)=U(±) [s(^'±)]"^ (28) 



and where U^^'^^ is the unitary matrix that diagonahzes [8^^'='=)]^^ (the overlap matrix, or 
its inverse, in the x-basis) 

U(±)t [s(^'±)]^' UW = [s(^'±)]^' (29) 

such that s(^'='=) is the diagonal matrix of eigenvalues and [s*^^'^^] ^ is the diagonal matrix 
of one over the square root of the eigenvalues. Note that in Eq. (29), the ±1 is independent 
of the ± which appears in the superscript (x^i)- Thus, 



S(<^,±) = 1 = X(±)^S(^'±)X(±); (30a) 
X(±)tx(±) ^ X(±)X(±)^ = X(±)X(±) = [S(^'±)] -\ (306) 

It can easily be seen that 

X(±)"' = U(±) [s(^'±)]^ U(±)^ (31a) 
XW =U(±)* [s(^'±)]^ UW. (316) 

II.B-2 Unitary Rotation Transformation 

In this section, it is assumed that the EOM in the 0-sct has been solved over the 
interval At so that the L>-matrix in the (p-set can be calculated at the R-end of the time 
interval. Thus, this section is concerned with the following derivative overlap matrices 
calculated at the R-end of the time interval: 
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~(±) r(±) 

D(^'±) =< \ (j) >; 

~(±) t(±) 



(32a) 
(326) 



The objective is to find a unitary transformation (at the R-end of the time interval) of the 
0-set to the ip-set, such that D^'^'='=^ = 1. Thus the transformation is represented by 



(33) 



or in matrix notation, this is described by Eq. (7) above. Note that the overlaps, given by 
this transformation, are related by 



(34) 



Clearly, if the 0-set is orthonormal, and if b is unitary, then the ip-set is orthonormal. 
Now, if Eq. (33) is used in Eq. (32), then 



j)i<p,±)^^{±) ^ g(^,±)^(±) 



(35) 



If we demand D(^'='=) = 0, and assuming orthonormality of the ip-set, then 



b(±) + ]3(<^,±)b(±) = (36) 

Now D*^"*^'^^ is anti-hermitian since the '0-set is orthonormal, as can be seen from Eq. 
(34), given the orthonormality of the 0-set. It is then known [7] that use of the Cayley 
decomposition to propagate b^ in time, will preserve unitarity. Therefore, if Eq. (36) is 
solved, using D*^'^'='=) at the R-end, then b^^^ will be calculated at the R-end such that 
D{^,±) ^ at the R-end. 



II. C EOM Invariance 



The objective of this section is to show the invariance of the EOM with respect to 
the orthonormalization transformation (X) and the rotation matrix (b) . This will be done 
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by first surmising a compact form for the full two-electron wave function (see Eq. (2)), 
and then deriving the EOM without assuming orthonormality or null D-matrices, for each 
of the three bases i/;). It will then be shown that the resultant EOM are invariant 

in form with respect to the linear transformations given by X and b. Note that, from 
hereon, unless otherwise noted, the ± notation will be dropped with the understanding 
that the distinction still applies. Collecting some pertinent formulas for reference, we show 
the following, which are inferred from Eqs. (4, 20, and 24) above: 



X - set 

*^(l,2) = ^(l)a(^)x(2); (37a) 

^<~X\H\~X> . (37c) 
— set 

vE'4l,2)=(^(l)a(^)0(2); (38a) 

#(^) = a(^)0; (386) 

h^"^) =< (^1^1 0> . (38c) 
•0 — set 

*^(l,2) = V^(l)aW 7/^(2); (39a) 

= ^{^)^. (395) 

hW=<^|^|^>. (39c) 



The actual orbital transformations have been given by Eqs. (7, 27). It must be shown that 
the EOM are invariant with respect to the two transformations. To this end, the EOM 
are derived below without assuming orthonormality or a null D-matrix. 

Begin by substituting Eq. (37a) into Eq. (2), multiplying from the left by x(l) and 
integrating over the coordinates of r(l)-the result is 

D Wa^x)^ + S^^^^a^^^X + S^^^a^^^^;^ = h^^^a^^^X (40) 
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Now, define a transformation of orbitals 



V' = x c 

(41) 

where 



C(±)=x(±)b(±); (42a) 

C(±)*=x(±)*b(±)*; (426) 

C(±)' =b(±)'x(±)'; (42c) 

cW=b(±)'xW"^ (42d) 

C(±)-^=xWb(±)\ (42e) 

Thus the transformation defined by C combines the orthonormahzation transformation 
with the rotation transformation. In this manner, the 0-set of orbitals may be dispensed 
with. Inverting Eqs. (41) results in 

(43) 

Then taking the time derivative of Eqs. (43), gives 



(44) 

X = C ^f+C-'^f 



Substituting Eqs. (41), (43), and (44) into Eq. (40) gives 



' <iP\ [if C"^ + if C"^] > a^^) C"^ ' < if\if > C-^a^^) C"^ if+ ^^^^ 

C^~' < if\if > C-^a(>^) [C V+ C- V ] = -i < ^\H\^ > C-^a^^) C"^ if. 
Now define 



i(^) = C-^a(>^)C-^ (46) 
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and then 

^C-ia(x)c-i + C-^a(^)C-^ + C-^a(^)6 \ (47) 

Also, it is clear that 

h(V') = h(^) C (48) 
Then multiply Eq. (45) from the left by C''' and use Eqs. (46) and (47) to produce 

DWaW^/;+S(^^a(^)V' + S(^)a(^)V' = h('^^a(^¥ (49) 

Note that 

SW=<V7|V'>; (50a) 

D(^)=<V^|V^> (506). 

A comparison of Eq. (49) and (40) demonstrates the EOM invariance. 
Now, using Eqs. (41) and (44), it is easy to see that 

D(^) = [dW C + C] (51) 
Requiring that D^'^) = results in an equation which defines C 

C + r(^) C = (52) 

where 

rW = sWd(>^) (53) 
In principle, r'^-''^) is antihermitian since S^-*^^ is hermitian and D^-^^^ is antihermitian. 
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III. USE OF THE ISOP ALGORITHM 



The ISOP algorithm[4] has been appHed to several problems [4,5] including the inte- 
gration of two coupled equations resulting from a time-dependent Hartree ansatz applied to 
iy2-[5] The ISOP algorithm is a two-time-step algorithm and involves (as indicated above) 
a time interval At, the retarded time (R), at the R-end of At, and the advanced time at the 
A-end. The present equations are nonlinear equations-this presents no essential difficulty 
however. It is possible, as we demonstrate below, to decompose the solution into steps so 
as to in effect render the equations linear from a computational point of view. The basic 
computational tool used in taking time derivatives in the ISOP is the Cayley formula. [7] 
The superscripts labeling ± and the orbital types are dropped in the following. 

First, rewrite Eq. (14) in matrix form as 



A = -ml (54) 

Then, the analysis easily gives 



Aa ~ [1+'- dt Ur] dt Ur] Ar (55) 

where the A-subscript stands for advanced (in time) and the /^-subscript stands for re- 
tarded (in time), as has already been used above. Also, 

Hk,j=«^k\H\^j» (56) 

and then 

[AA-AR]/At. (57) 
The time derivative of the orbitals is evaluated in this manner also; e.g. for the X"Set 
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(58) 



The time advance of the matrix C may be evaluated (see Eq. (52)) in a similar manner 
(when needed) via 



Ca ^ [1 + ^ dt Tr] '[1 - i dt Tr] Cr (59) 

It is useful to express the ilf-matrix elements (Eq. (55)) in terms of the /i- matrix 
elements. To do so, define 



hp{qr)s =< V'pI VI V's > (60) 



Then, (with Hjk = Hkj), 



Hii = ; 

Hi2 — [/ii(ii)2 + ^i(i2)i] ; 
Hi3 = /ii(i2)2 ; 

-H^22 = - [/i2(ll)2 + 2 /l2(12)l + /il(22)l] ; (61) 

H23 = [^1(22)2 + ^2(12)2] ; 

-f^33 = /i2(22)2 ; 

Hii = - [/i2(ll)2 - 2 /i2(12)l + ^1(22)1] • 

In order to use the ISOP to integrate the orbital equations, the singlet and triplet 
cases are considered separately again. The ISOP may now be applied directly to Eq. (19). 
Using the definitions of To{j) and Vb(j) immediately below Eq. (3), Eq. (19) may be 
written as (the explicit reference to particle coordinates is dropped) 

* = -z To* - z M * (62) 

where 

M = IVb + ho + V (63) 

such that ho is the same as h in Eq. (23), except with H{1,2) replaced by Ho{l,2) (see 
Eq. (3)). Also V is the same as h except V{1, 2) replaces H{1, 2). 
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Then, finally, application of the ISOP results in 



1 - ^idt Tq. . . 1 - hdt To ^ 

where 

Mi,ot = +lidt M 

1 . (65) 
Mtop = --idt M 

Assuming that the X"Set is orthonormal and that D^^) = 0, as would be the case for 
a completely accurate time propagation and spatial derivative and integral evaluation, the 
basic algorithm proceeds as follows (reintroduce the orbital labels): 

Basic Algorithm 

Step 1: specify XR{t = 0),A'i\t = 0) = 1; 
Step 2: calculate aS^^ using Eqs. (22,23); 
Step 3: calculate h^"* using Eq. (37c); 
Step 4: calculate H^'* using Eqs. (59); 
Step 5: calculate yl^^ using Eq. (55); 
Step 6: calculate using Eq. (39b); 
Step 7: calculate using Eq. (62); 
Step 8: calculate a^'' using Eqs. (22,23); 

Step 9: calculate xa using Eqs. (37b) after having taken the inverse of a^^; 
Step 10: calculate x using Eqs. (58); 
Step 11: calculate S^^ using Eq. (25); 
Step 12: calculate D*^^) =< x\ X >'i 

Step 13: if S^"* ft; 1 and if D*^^) ^ 0, continue, otherwise, jump out of this algorithm; 
Step 14: let the advanced time become the retarded time for another time interval 
and then go to step 3 and continue. 



If an exit of the basic algorithm occurs because of the lack of orbital orthonormality or 
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the lack of a null D-matrix, the "Correction Algorithm" is implemented, for the particular 
At in question, as follows: 

Correction Algorithm 

Step 1: assume = 1; 
Step 2: calculate Tji using Eq. (53); 
Step 3: calculate using Eq. (59); 
Step 4: calculate i/ja using Eqs. (41); 
Step 5: calculate using Eq. (46); 
Step 6: calculate using Eq. (48); 

Step 7: let the advanced time become the retarded time for another time interval and 
then go to step 4 of the "Basic Algorithm" and continue (letting i/j —>■ x)- 
Thus, the solution may be efficiently propagated in time in a way that guarantees orbital 
orthonormality and a null D-matrix (non-reduntant time evolution). 

IV. CONCLUSIONS 

Now returning to the TDHF[9], the one-electron orbitals are assumed orthonormal and 
the derivative matrix elements are set to zero. This is possible in the context of the TDHF 
because the EOM are derived using the Dirac-Frenkel variational principal (DFVP). [10,11] 
In fact, orthonormality of the one-electron orbitals and zero D-matrices are not restrictive 
in the context of the TDHF. This can be demonstrated by a set of guagc transformations 
in the manner of Refs. (5) and (9). The key idea of the TDHF derivation using the DFVP 
is the independent variation of the one-electron orbitals. [9] The authors of Refs. (2) and 
(3) also use the DFVP to derive their EOM for a Hartree-like ansatz which is used for 
distinguishable particles. The present derivation of a set of exchange equations proceeds in 
the manner of Ref. (5)-that is by projections. It is necessary to have a way of determining 
the D- matrices. This is so because the time derivative that appears inside the spatial 
integral is at the same time-step as the time-derivative of the Schrodinger equation itself. 
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This results in the necessity of a self-consistent procedure at each time-step-an operation 
that needs to be avoided if possible. If projections are used to derive the two coupled 
equations for the one-electron orbitals, the D-niatrices can not be set to zero without a 
"hidden symmetry" analysis similar to that given above. This in fact has been done. [12] 
The equations that result from that procedure are distinctly different from the TDHF 
equations. 

The present work is an alternative derivation of a set of time-dependent exchange 
equations (TDEEs), which follows immediately from the MCTDH theory of [2,3]. The 
present equations are derived by a projection procedure, but are exactly the same equa- 
tions that would be produced by the DFVP. This all results from the reduntancy which 
follows from the use of the purely time-dependent coefficients multiplying each two-electron 
configuration. 

An algorithm for the solution of the new exchange equations has been given using the 
ISOP method and a procedure was described to enforce orbital orthonormality and a null 
D-matrix. This allows for a non-reduntant time evolution. 



ACKNOWLEDGEMENTS 



The research was supported by the Army High Performance Computing Research Cen- 
ter and the US Army, Army Research Laboratory (DA AH04-95- 2-0003/ contract number 
DAAH04-95-C-0008), by NSF CREST grant HRD-9707076, and by the Lawrence Liver- 
more National Laboratory Research Collaboration Program for Historically Black Colleges 
and Universities and Minority Institutions. The author would like to acknowledge useful 
conversations with H.-D. Meyer, Burke Ritchie, Merle Riley and Mario Encinosa. 



19 



REFERENCES 



[1] C.F. Barnett, "Atomic Collision Properties" in A Physicist' s Desk Reference : 

The Second Edition of Physics Vade Mecum, American Institute of Physics, New 
York, N.Y., ed. H.L. Anderson (1989), p. 92. 

[2] H.-D. Meyer, U. Manthe, and L.S. Cederbaum, Chem. Phys. Lett. 165, 73 (1990). 

[3] U. Manthe, H.-D. Meyer, and L.S. Cederbaum, J. Chem. Phys. 97, 3199 (1992). 

[4] B. Ritchie and M.E. Riley, Sandia Report Sand97-1205, UC-401 (1997). 

[5] B. Ritchie, C.A. Weatherford, International J. Quant. Chem. S70, 627 (1998). 

[6] E.K.U. Gross, E. Runge, and O. Heinonen, Many — Particle Theory, Adam Hilger, 
New York (1991). 

[7] M.D. Feit, J. A. Fleck, and A. Steiger, J. Comput. Phys. 47, 412 (1982). 

[8] A. Szabo and N.S. Ostlund, Modern Quantum Chemistry, McGraw-Hill, New York 

(1989), pgs. 142-145. 
[9] A.K. Kerman and S.E. Koonin, Ann. Phys. (N.Y.) 100, 332 (1976). 
[10] P.A.M. Dirac, Proc. Cambridge Philos. Soc. 26, 376 (1930). 

[11] J. Frenkel, Wave Mechanics : Advanced General Theory, Clarendon (Oxford) (1934). 
[12] M. Riley, B. Ritchie, and C.A. Weatherford, unpublished work. 



20 



